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ABSTRACT 


The experimental verification of active control methods for vibration suppression 
of large flexible structures in space is essential for precision optical and military payload 
operations. The Flexible-Spacecraft Simulator (FSS) at the Naval Postgraduate School is 
designed for testing such control designs. The experimental setup simulates the pitch 
axis motion of a rigid body spacecraft with a flexible antenna support structure connected 
to a rigid reflector. A twenty-four state finite element analytical model is used to 
characterize the flexible appendage. Piezoelectric sensors and actuators are used for 
feedback control for vibration suppression. In addition, an external infrared camera 
provides direct feedback of the flexible structure’s elbow and tip displacements and 
rotations. A Multiple-Input-Multiple-Output (MIMO) linear quadratic gaussian (LQG) 
controller is designed using linear quadratic regulator (LQR) optimal control theory and 
an optimal Kalman estimator as the state observer to meet desired performance 
specifications. The objective is to minimize the motion of the reflector. 
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[. INTRODUCTION 


The rapid advancement of space sensor technology and the requirement for 
nanometer accuracy has created a need for minimum vibration noise, high precision 
mounting platforms from space. The consumer’s demand for high resolution imagery, 
large space structures, and precise position and attitude determination has produced a 
large amount of advanced research in the field of structural control. Concurrent to the 
rapid increase in accuracy and resolution of sensor payloads is the need to isolate that 
payload from the disturbances inherent to any spacecraft. 

For single event upsets, the control problem lies in damping out vibrations 
incurred through the disturbance. Slew maneuvers, thruster firings, micro-meteorite 
collisions are all examples of SEU disturbances. Multiple event upsets, periodic or 
continuous disturbances, come in the form of high frequency vibrations from an 
internally mounted momentum wheel assembly, a cryogenic cooler apparatus, or any 
vibrating machinery that is part of the main spacecraft bus. The method of control for 
this type of disturbance is vibration isolation, and its control synthesis is inherently 
different from the damping control problem. In order for a precision sensor platform to 
be stable in the space environment and exhibit robustness in disturbance rejection, both 
types of control must be combined and used. 

Optimal control is an excellent method to achieve both control specifications. 
For these types of applications, smart structures are a promising technology. In general, 
smart structures are the system elements that sense the dynamic state and change the 





system’s structural properties, such as its natural frequencies and its damping, to meet 
given performance objectives. 

There are several types of embedded sensors and actuators which can be used for 
vibration and structural control. The embedded sensors are piezoelectric deformation 
sensors, strain gages, and fiber optic sensors. The embedded actuators are piezoceramic 
wafers, electrostrictive ceramic wafers, piezoceramic polymer film and shape memory 
metal wires. Piezoceramic sensors have a high strain sensitivity, a low noise baseline, low 
to moderate temperature sensitivity, and an ease of implementation. Piezoceramic 
actuators have high stiffness, sufficient stress to control vibration, good linearity, 
temperature insensitivity, are easy to implement, and minimize power consumption. 

Conventional control methods have worked well in the past, but new design 
methods are required to obtain improved performance and robustness characteristics from 
the structural control system in order to satisfy future design specifications. Positive 
position feedback (PPF) and velocity feedback are two proven methods of structural 
control that work well with piezoceramic actuators and sensors. With a multiple-input- 
multiple-output (MLMO) control system, linear quadratic control methods are a preferred 
choice. 

Linear quadratic control is based on full state feedback. If all the states are not 
known, an observer is inserted into the loop to estimate the unknown states. Linear 
quadratic controllers using state estimate feedback are optimal for the nominal plant 
models but the performance may be far from satisfactory in real life due to plant 
uncertainties, unxnodeled plant dynamics, and sensor noise that is unaccounted for in the 


system compensator. Attractive passband robustness properties of full state feedback 



optimal quadratic designs may disappear with the introduction of a state estimator. The 
Linear-Quadratic-Gaussian methodology provides an mtegrated frequency domain and 
state space approach for design of IVOMO control systems. The advantages of the 
methodology lie in its ability to directly address design issues such as stability robustness 
and evaluate the trade-off between performance and allowable control authority. 
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II. SCOPE OF THESIS 


The objective of this research is to achieve vibration suppression of a flexible 
antenna structure and also isolate vibration from the spacecraft main body. The flexible 
spacecraft simulator (FSS) in the Space Systems Dynamics and Control laboratory will be 
used to experimentally verify the data obtained from this analysis. The inherent 
spacecraft noise is in the form of vibrations from a continuous source, a momentum 
wheel mounted inside the spacecraft. The use of piezoceramic sensors and actuators to 
negate the disturbance vibration will be used. In addition, an external infrared optical 
sensor will be employed to provide structure position and velocity information to the 
optimal controller. The control technique we implement will be Optimal Linear 
Quadratic Gaussian (LQG) controller. 

The approach will be to first characterize the flexible beam structure using finite 
element analysis and verify the results with elementary modal analysis. Using the 
modeshape information from this analytical model, optimum placement of the 
sensor/actuator pairs will be determined. Next is to construct an optimal control system 
using Linear Quadratic Regulator (LQR) optimal control with multiple sensors and 
actuators. Once the LQR control gains are satisfactory, an optimal observer will be added 
to the feedback path for smoothing and to estimate the unknown states plus filter out any 
high frequency noise. 

The use of multiple sensors and actuators will enable the controller to quickly 
dampen out higher order modes of the reoccurring disturbance. Using SIMULINK, a 
vibration isolation system simulation will be designed, implemented and tested. The 



optimal controller will be compared against known effective control schemes such as 
Positive Position Feedback (PPF) and derivative feedback control. 

The follow on research will be to implement the LQG controller on the real time 
controller in the spacecraft dynamics and control laboratory, AC-100 from Integrated 
Systems, Inc. 
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III. CONTROL THEORY 


A. BACKGROUND 

State-space methods are the cornerstone of modem control theory. The essential 
difference between modem control theory and classical control theory is the 
characterization of a process by differential equations instead of transfer functions. In the 
modem approach, processes are represented by systems of coupled, first-order differential 
equations. In principle, there is no limit to the order of the system and in practice the 
limit is the availability of the computer software to perform the required calculations 
reliably. 

B. STATE SPACE REPRESENTATION 

Differential equations can be expressed as a set of simultaneous first-order 
differential equations. They are represented in state-variable form as the vector equations 


X 

y 


f(x,u) 

h(x,u) 


(3.1) 


where the input is u, and the output is y. 

The column vector x is called the state of the system and contains n elements for a 
nth-order system. For mechanical systems, the states usually consist of the positions and 
velocities of the separate bodies. The vector function f relates the first derivative of the 
state to the state itself and the input u. Parameters in the function f could be dependent on 
time as well as position and velocity, but by and large, for structural systems the 




parameters can be considered time-invariant. Linear approximations will be used for 
nonlinear systems in design and analysis. For the linear case, equation (3.1) can be 
written as 


X = Ax + Bu 
y = Cx + Du 


(3.2a,b) 


where A is an nxn system matrix, B is an nxm control matrix, C is an Ixn observer matrix 
and D is an Ixm feed-through matrix(direct transmission matrix). Under most 
circumstances D is normally the zero matrix. 

C. TRANSFORMATION OF STATE 

For a system described by equations (3.2), that description is not a unique state 
description of the system. Matrix algebra enables a linear transformation of state without 
changing the basic properties of a matrix. Consider a state vector z where 

x = Tz (3.3) 


and T is a nonsingular matrix. Substituting equation (3.3) into equation (3.2a), a linear 
transformation of state is performed 



x = Tz = ATz + Bu 


z = T ‘ATz + T 'Bu 
z = AtZ+ ByU 

where (^.4) 


At = T ‘AT 


Bt = T 'B 


Substituting equation (3,3) into equation (3.2b) 


y = CTz+ Du 


= CtZ + Du 




(3.5) 


It is sometimes convenient to transform a physical system model into its modal 
canonical form, also known as the Jordan canonical form. This state description 
decouples the coupled first order equations into n independent first-order equations, 
providing that the system matrix is diagonalizable. This description has many advantages 
such as the degree of control authority the input has on each mode, the observability of 
each mode, and the damping ratio of each mode provided it is a damped system. 

D. OPTIMAL CONTROL THEORY 

An effective and widely used control technique of linear control systems is the 
optimal Linear Quadratic Regulator (LQR). Provided the full state vector is observable. 
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this method can be employed to fit specific design and performance criteria. A quadratic 
cost function, based on the Bolza problem, is used to minimize the performance index, J. 
The general form for the LQR is 

J = J(x^Qx+u'Ru)r/t (3.6) 

where Q is the state weighting matrix and R is the control weighting matrix. The 
necessary conditions for the optimal solution are that Q must be symmetric and positive 
semi-definite and R must be symmetric and positive definite. The solution to the LQR 
problem results in the optimum full state feedback gain matrix 

K = R-'B’M (3.7) 

where M satsfies the matrix Riccati equation 

-M = MA + A'M-MBR“'B'M + Q (3.8a) 

For a time-infinite solution, M is set to zero. The control input is then 

u = -Kx (3.8b) 

If the full state vector is not available, in order to use LQR control theory, it is 

necessary to design a state estimator or observer. The general form for an observer is 


given by 



X = Ax + Bu + L(j/ - Cx) 


(3.9) 


where x is the estimated state vector. The inputs to the observer are the outputs from the 
plant and the control inputs to the plant. 

By solving for the optimum observer gain, L, the observer is known as a Kalman 
filter. The general form of the dynamical equations used in the Kalman filter synthesis 


X = Ax + Bu + Fv 
y = Cx + Du + w 


(3.10) 


where F is the plant uncertainty matrix, w is the state noise vector, and v is the sensor 
noise vector. Both v and w are considered to be white noise processes. White noise 
processes have the property of having a mean value of zero over time. This assumption 
does not always hold true in the real world but it simplifies the analysis considerably. 
The solution for the optimal observer gain is given by probability theory. The Kalman 


observer gain, L, is given by 




L = PC'W ' 


(3.11) 


where P is taken from 

P = AP + PA'-PC'W^'CP + FVF' (3.12) 

The process noise covariance matrices V and W are given by 

£{v(0v'(0} = V(0<5(t-r) 

£{v(/)w'(0} = X(/)<y(/-r) (3.13) 

£M/)w'(0} = 'V(05{t-r) 

X(t) is the system cross-covariance matrix and is a function of the correlation of sensor 
noise to plant noise and under most circumstances it is normally zero. The symbol E{ } 
denotes mathematical expectation, in other words the average computed in the 
probabilistic sense. 

E. CONTROLLABILITY AND OBSERVABILITY 

Controllability and its dual, observability, are system parameters that indicate the 
extent that the states are controllable via the control matrix and for observability, the 
number of states that are observable via the sensors. They are two mutually independent 
functions of the system plant matrix, control matrix and the observer matrix. One simple 
test for controllability and observability is to form their respective matrices, the 
controllability matrix (Me) and the observability matrix (Mo). Provided each matrix is of 
full rank, that indicates that the states are independent of one another and can be 
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controlled or obsei-ved. The degree of controllability or observability cannot be 
determined from this test; however converting the system to Jordan canonical fonn will 
indicate the relative control authority to the individual modal states and also the relative 
observability for each modal state. 

The following are the equations for the respective matrices: 


M, =[b AB A"B A'-'b] 

M,=[C’' A''C^ (A^)"C'' ••• (A^)"-'C’'f 


(3.14) 


A complete treatment is given in several textbooks, such as [Refs. 1,2]. These 
concepts will play a fundamental role in the design of the control system for the flexible 
structure. The whole concept of control design is a compromise between accuracy, 
complexity, and robustness. In particular, as this thesis will show, a very accurate model 
of the system might be not only too complex, but also it might be bound to be 
uncontrollable or unobservable, fhis is because an exact model is likely to contain a 
number of states, i.e., modes, which are hard to completely control with the given input 
signals. Likewise these same states might not be easily observable. When this is the 
case, the designer would look for a simpler model, where all the states are controllable 
and observable, and (hopefully) the effect of the neglected states will be within the noise 
level. For this problem the two matrices discussed above (controllability, observability) 
give a clear indication on the characteristics of the model. 
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IV. PIEZOELECTRIC THEORY 


A. HISTORY OF PIEZOELECTRICITY 

Piezoelectricity is a phenomenon that describes certain materials that generate 
electricity when a mechanical stress is applied, known as the direct effect, and also, that 
generate a mechanical stress when an electric field is applied on them. Piezoelectricity 
occurs naturally in some crystalline materials and can be induced in other polycrystalline 
materials through a process known as “poling”. The material’s crystal lattice structure 
may be poled by the application of a large electric field, usually at high temperature (see 
Fig. 4.1). [Ref. 3] 



5^ 
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\ / / 
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Figure 4.1. [Ref. 3) Crystalline Lattice Structures 
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Subsequent application of an electric field will produce additive strains locally 
which translate into a global strain for the material. The direct piezoelectric effect has 
been used for a long time in sensors such as accelerometers. Use of the converse effect, 
however, until recently has been restricted to ultrasonic transducers. Barium titanate, 
discovered in the 1940s, was the first widely used piezoceramic. Lead zirconate titanate 
(PZT), discovered in 1954 [Ref 4], has now largely superseded barium titanate because 
of its stronger piezoelectric effects. Only recently have researchers in the area of 
structural control taken notice of the very desirable features of piezoelectric actuators and 
have started using them for many structural control applications. Piezoceramics are 
compact, have a very good fi-equency response and can be easily incorporated into 
structural systems. Actuation strains on the order of 1000 pstrain have been reported for 
PZT material. Strains are non-dimensional numbers which relate the change in length to 
the original length for a given impetus. Within the linear range they produce strains that 
are proportional to the applied electric field/voltage. These features make them very 
attractive for structural control applications. 

There are several methods to model the constitutive behavior of piezoelectric 
materials. The most popular is the macromechanical approach, it provides the 
relationship between the electrical and mechanical effects in a manner that can be 
incorporated onto typical isotropic or orthotropic structural materials. For linear 
piezoelectric materials, the interaction between the electrical and mechanical variables 
can be described by linear relations of the form 
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S, =s^T,+d„,E^ 
D^=d^J,+£l,E, 


(4.1) 


The mechanical variables are the stress T and the strain S, and the electrical variables are 
the electric field E and the eleetric displacement D; s is the compliance, d is the 
piezoelectric constant and e is the pennittivity. The first equation describes the converse 
piezoelectric effect, and the second equation describes the direct effect. ITie stress and 
strain are second order tensors, while the electric field and electric displacement are first 
order. 

Figure 4.2 shows the typical coordinate system used to represent a poled 
piezoelectric. The 3-axis is in the direction of the initial polarization. 


■7^ 



indicates poling direction 


Figure 4.2. [Ref. 3] Typical Piezoceramic Coordinate System 


The 1- and 2-axes are arbitrary in the plane perpendicular to the poling direction. The 1 - 
and 2-axes are arbitrary because a poled piezoelectric is transversely isotropic in the 1-2 
plane. The equations above written explicitly in matrix form are 
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0 0 0 0 0 



0 

0 

0 

0 





(4.2) 


Where S| through S3 are the normal strains, S4 through 85 are the shear strains, T| 
through T3 are the normal stresses, T4 through Te are the shear stresses, D| through D3 
are the electric displacements and E| through E3 are the electric fields associated with the 
given coordinate system. 

The piezoelectric constants that are of most interest from a structural standpoint 
are the d constants. These constants relate the strain developed in the material to the 
applied electric field; obviously, the highest value of these constants is desirable. The 
constant relates the strain in the 3-direction to the field in the 3-direction.; similarly, the 
dsi relates the strain in the 1-direction to the electric field in the 3-direction. The electric 
field is voltage applied across the piezoelectric divided by its thickness. It is important to 
point out that rfjj is usually a positive number and d^i is a negative number. This means 
that a positive field (i.e., a field applied in the poling direction) will produce a positive 
mechanical strain in the 3-direction and a negative strain in the 1-direction. 

There are many applications for the use of PZT actuators and they can be divided 
into two general categories: the linear actuators and the actuators used for structural 
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control applications. In the first category, the PZT actuators, arranged in the form of 
stacks, are used in a fashion similar to shakers or conventional hydraulic or electrical 
actuators. Due to the fact that the stroke is severely limited, in static applications they 
are used only for micropositioning work. In structural control applications, the actuators 
are typically embedded or bonded to the surfaces of the structure and they apply localized 
strains that can be used directly to control structural deformations. A typical 
arrangement is shown in Figure 4.3; the two actuators are on the upper and lower surfaces 
of the structure and are actuated out of phase (the upper expands and the lower contracts), 
which creates a moment on the structure. [Ref 3] 



In the remainder of this chapter we will apply piezoceramic actuators and sensors 
structiual control elements. 
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B. FLEXIBLE SPACECRAFT SIMULATOR PIEZOCERAMIC 
ACTUATORS AND SENSORS. 

Piezoceratnic sensors and actuators are located on the flexible appendage as 
shown in Figure 4.4. The piezoceramic wafers are bonded to the surface of the flexible 
arm and the voltage developed from the sensors is fed to the actuators by way of the 
designed control system. 



Figure 4.4. Piezoceramic Actuators and Sensors Mounted on the FSS Flexible Beam 


Figure 4.5 illustrates the orientation of a piezoceramic wafer on an arm and the alignment 
of its axis that describes the electro-mechanical relationships. 



20 





The piezoceramic wafers in a sensory mode produce a charge between their 
electrodes that is directly proportional to the lateral strains. It is given by 


Q=AEd^^{e, 


where A is the lateral area of the piezoceramic water, E is Young’s modulus of the wafer, 
dj! is the lateral charge coefficient, and ei and £2 are the strain values in the lateral 
directions respectively. The capacitance for a piezoceramic wafer as shown in Figure 4.5 
is given by 


C = 


DA 


(4.4) 


where D is the dielectric constant of the piezoceramic and t is the thickness of the wafer. 
The voltage V produced by a sensor under strain is given by 


V 


(4.5) 


When using piezoceramic wafers as actuators, the attachment geometry is similar to the 
sensor geometry shown in Figure 4.5. The control voltage , e^, is applied to the wafers 
and the lateral strain that is developed can act to control the bending of the beam. The 
electric field that is developed by the wafer is given by 
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(4.6) 


Care must be taken not to induce a strong electric field that is opposed to the 
piezoceramic’s poling direction as that can damage the material by depolarizing it. 
Typical field limits by most materials are between 500 and 1000 volts/mm. 

Table 4.1 is a good reference for the material properties of piezoceramics. For 
this thesis. Navy Type 11 PZT piezoceramic element is used, which is material EC-66 in 
Table 4.1. It is a above average piezoceramic material which boasts an excellent Lateral 
Strain coefficient (t/31) and exhibits good linearity properties between applied voltages of 
±150 volts. 
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V. FINITE ELEMENT ANALYSIS 


A. BACKGROUND 

The need for characterization of complex structures and the advent of the modem 
day computer has given rise to a method of analysis known as the finite element method. 
The idea behind the finite element method is to provide a formulation which can exploit 
digital computer automation for the analysis of irregular systems. To this end the method 
regards a complex structure as an assemblage of finite elements, which every such 
element is part of a continuous structural member. By requiring that the displacements be 
compatible and the internal forces in balance at certain points shared by several elements, 
where the points are known as nodes, the entire stmcture is compelled to act as one entity. 

The finite element method of analysis, while considering a continuous structure, is 
a discretization problem in that it expresses the displacements of the continuous structure 
in terms of a finite number of displacements at the nodal points multiplied by the 
interpolation functions. The advantage of the finite element method over any other 
method is that the equations for the system can be derived by first deriving the equations 
for a single element and then assembling the individual elements’ equations by using 
constraint conditions. The displacement at any point inside the element is obtained by 
means of interpolation, where the interpolation functions are generally low-degree 
polynomials and they are the same for every element. [Ref 5] 
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ELEMENT STIFFNESS MATRIX 


Meirovitch [Ref. 5, pp. 303-304] uses the direct method for derivation of the 
elemental stiffness matrix. The stiffiiess matrix relates a displacement vector to a force 
vector. Using the element in Figure 5.1, for uniform bending stiffiiess, the differential 
equation for the displacement w(xj is 


d*w{x) 


(5.1) 


where E is Young’s modulus for the material and / is the area moment of inertia for the 
beam cross-section. 
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By integrating four times, we have the elemental displacement at any point 


w(x) = -^c,x^ +c^x + c^ (5.2) 

with C] through C4 as the constants of integration, determined from the boundary 
conditions 


w(0) = w, 


dw{x) 

dx 


w{h) = w, 


dw{x) 

dx 


(5.3) 


w\ and W2 are the nodal displacements and 0] and di are the nodal rotations at the 
endpoints of the element. By solving equation (5.2) into (5.3) the corresponding 
constants of integration are solved as 


•I “ -2 vi'2 +/i6^) Cj = -^(-3^1 - 2h0, +3wj -/i^j) 


(5.4) 


Hence, introducing equation (5.4) into (5.2), the expression for the bending displacement 
is determined to be 
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*(i 


The bending displacement is related to the inter-element nodal forces/1 ^ y 3 , and/4 as follows 



(5.6) 


Combining equations (5.6) and (5.5) yields 


/, =p-(12H', +6/20, -12^2 +6/?(9j) 
fj = ■^(6>»'| +4/20, -6^2 + 2/202 ) 

/j =-^(-12 h', -6/20, +12w2 -6/202) 

El 

/a =-^(6w, +2/20, -6^2 +4/202) 


This can be written in matrix form as 


wW={/} 


(5.8) 
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where the nodal displacement vector and the nodal force vector are given by 



w, 



M = 


{/}=• 

f2 

\A 




u. 


and the element stiffiiess matrix is given by 


12 6 h -12 6 h' 

£/ 6/t 4h^ -6h 2h^ 

-12 -6h 12 -6h 


6 h 2h^ -6h 4h^ 


(5.10) 


C. ELEMENT MASS MATRIX 

Equation (5.5) can be written as a function of the nodal displacements and 
interpolation functions, L(x) 


r (xV fxV] 

i3(^)=| 

\(xY fxVl 



\ h) \ hJ j 

\hJ KhJ J 




\_h ^\h) UJ 

•^4 (^) = 

UJ UJ J 
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Equations (5.11) and (5.12) can be expressed in matrix form 

w(x,0 = {i(x)}^{x(0} (5.13) 

with {LfxJ} as a four-dimensional vector of the interpolation functions and {x(t)} a four¬ 
dimensional vector of nodal displacements. The element kinetic energy has the form 

W = ^ dx^\{w{t)Y[m\{w{t)} (5.14) 

where the 4x4 mass matrix is given by 

[m]= |'/?j(x){L(x)}{L(x)}’^£)!x (5.15) 


Inserting the interpolation functions vector, equation (5.12), into equation (5.15) 
and integrating over the element length, the elemental mass matrix becomes 
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■ 156 22h 54 -Uh 

mh 22h 4h^ I3h -3h^ 

4^ 54 \3h 156 -22h 

-\3h -3h^ -22h Ah- 


(5.16) 


D. ASSEMBLING THE GLOBAL MASS AND STIFFNESS MATRICES 
INTO THE COMPLETE SYSTEM 

The next step is to assemble each element’s mass and stiffiiess matrices into the 
global system mass matrix and the system global stiffness matrix, respectively. The 
assembling procedure is an algebraic process of adding the overlapping element matrices 
together. Consider the 4x4 element matrix partitioned into a 2x2 matrix of 2x2 sub¬ 
matrices 






(5.17) 


The global system matrix is required to be symmetric and is constructed by adding the 
element matrices along the diagonal 




(5-18) 


For example, with three elements this matrix becomes 
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(5.19) 


0 0 ] 



The method is the same for both the stiffiiess and the mass matrices. In order to 
satisfy the boundary conditions for a fixed-free system, the displacement and rotation at 
the root of the flexible appendage must be zero. This condition is satisfied by eliminating 
the first and second rows and the first and second columns from both the global mass and 
stiffriess matrices, respectively. 

E. PIEZOELECTRIC FINITE ELEMENT MODELING 

The piezoceramic elements used in the structure must be incorporated into the 
finite element model. The basic equations for both piezoceramic actuators and sensors 
are the same as for ordinary structural elements discussed in the previous section. There 
is a need to compensate for the piezoceramic displacement from the center of the beam 
and Figure 5.2 illustrates this concept. In addition, the electro-mechanical relationships 
must be taken into account for implementation into an analytical model suitable for 
control design. 



Figure 5.2. Piezoceramic Element Placement on Beam 
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From equation (4.1) the general relationship for the electro-mechanical coupling is 


given by 


Using the fact that the elastic constant for piezoceramic material, s, is the inverse of its 
Young’s modulus, Ep, this equation can be written as 


-d,,E^ E^ 


(5.21) 


I'he next step is to set up the equation for the elemental potential energy, U. 


-1/= 1 J(-7;5, +D,E,)dV 


(5.22) 


where the two terms in the integral represent mechanical energy (T’lS'i) and electrical 
energy Using Wp as the width of the piezoceramic wafer, this equation can be 

rewritten as 



(5.23) 


33 




The strain, using small angle displacement theory, Si, can be written as 


S, =£, 


(5,24) 


substituting equation (5.21) into equation (5.23) 


'yiE,Yh;^d,‘E, d„E,]\E,\ 




= \[{el + 2d,,E^E,e,-E^e,'^]dxdz 


then, using equation (5.24) 


U = -w„ f \\(£l -d,,^EJE,^ +2d,,EE,z^~E, 
2 ^ i \ ^ ^ ^ 


Rewriting equation (5.13) in terms of a summation 


w(x,0 = Z3>,(:r)9,(0 


(5.27) 
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where O is the vector of interpolation functions or “modeshapes” and q is the nodal 
vector fiom equation (5.9a). Inserting this equation along with the interpolation 
functions, equation (5.12), the general form of the energy equation is 


-U = —ye^ -q^be-~q^k^q 


where 


w hj ^ 

= —(^3 e = t^E 

^ f/- W , 


.2 . f/lrrf'<I>,(x)rf^<t>,(A:) _ 

C+Q,+-r —71- 

3 ^ dx dx 


(5.28) 


Substituting the interpolation functions from equation (5.12) into the b vector, we get 


=0 



and into the piezoceramic elemental stiffness matrix, [kp] 
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(5.30) 


The piezoceramic elemental mass matrix is derived in the same fashion as an 
ordinary structural element with the exception that the piezoceramic material linear mass 
density be used instead of the ordinary material density of the structure to which it is 
bonded. 

For structural elements that have piezoceramic material bonded to them, their 
respective mass and stif&iess matrices are given by the simple addition of the beam 
elemental matrices and the piezoceramic elemental matrices which are then assembled 
into global mass and stif&tess matrices as in the previous section. 


=[ml 





( 5 . 31 ) 
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VI. ANALYTICAL MODEL 


A. FLEXIBLE SPACECRAFT SIMULATOR ARM 

The Spacecraft Dynamics and Control Laboratory Flexible Spacecraft Simulator 
(FSS) simulates attitude motion about the piteh axis of a spacecraft. As shown in Figure 
6.1 it consists of a single degree-of-ffeedom rigid central body, representing the 
spacecraft central body, and a multiple degree-of-freedom flexible appendage, 
representing an antenna reflector with a flexible support structure. 



Piezoceramic sensors and actuators are used to provide active damping to the 
flexible support structure. The entire system is floated on air pads over a finely ground 
granite table to simulate a microgravity environment. The central body has two thrusters 
and a momentum wheel as its actuators. The flexible appendage has two stacked 
piezoceramic pairs as actuators. The first pair is located at the base of the arm assembly 
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and the second is located at the base of the forearm near the structure’s elbow as in 
Figure 6.2. 



Figure 6.2. FSS Experimental Setup for LQGA^ision Server Testing 

A finite element model will be used to simulate the flexible appendage. The 
elemental mass and stiffiiess matrices are given by equations (5.16) and (5.10), 
respectively. Using equation (5.18) to assemble both the global mass matrix, M, and 
global stiffness matrix, K, and by the use of Newtonian physics, the equations of motion 
can be derived. Figure 6.3 reveals the finite element model used for this research. 



Figure 6.3. Finite Element Model of the FSS Flexible Appendage 
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Table 6.1 gives the material properties and parameters used on this model. 


7075 T-6 Aluminum Beam 

Quantity 

Description 

Units 

Value 

tb 

Beam thickness 

meters 

1.5875x10’ 

Wb 

Beam width 

meters 

0.0254 

Pb 

Beam density 

kg/m^ 

2800 

Eb j Young’s Modulus 

N/m^ 

1.029x10’ 

Navv Tvne II PZT 

Quantity 

Description 

Units 

Value 

dll 

Lateral strain coefficient 

m/V or Coul/N 

-1.8x10'“ 

Ep 

Young’s Modulus 

N/rn^ 

6.3x10'“ 

v 

Poisson’s ratio 

N/A 

0.35 

D 

Absolute permittivity 

Farad/m or N/V^ 

1.5x10“ 


Table 6.1. Model Material Properties 


B. EQUATIONS OF MOTION 

In this section we derive the equations of motion for a twelve degree-of-freedom 
dynamic model. The derivation goes along the lines of a two-degree-of-freedom system 
as shown in Figure 6.4. 
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From Newton’s Law 


= ma 


(6.1) 


applied to each mass we obtain 

/n,Xi =-i|X|+^2(>^2-^i) - x^) + F (6.2) 


Using a vector 


, and grouping like terms 



where M is the mass matrix, K is the stiffiiess matrix, and F is the force vector. 

This type of dynamic representation relates directly to finite element analysis with 
appropriate modifications to account for piezoceramic sensors and actuators in the finite 
element model. 

Using the principle of virtual work and using the Lagrangian method of analysis, 
the kinetic energy for a piezoceramic element can be written as 
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(6.4) 


where is the mass per imit length for the piezoceramic and vi' is the velocity at any 
point on the element. An alternative form of equation (6.4) is 


(6.5) 


where 


[m^]. = \p^<j>Xx)tj)j{x)dx 


(6.6) 


given (f) the interpolation functions from equation (5.12). The total kinetic energy for a 
beam mounted piezoceramic is then 


r = 


-q’^Mq 


(6.7) 


where the mass matrix is given by equation (5.31). 


The Lagrangian function, L, is given by 
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L = T-U 

= + --q^'Kq 


(6.8) 


which is used in the equations of motion 


(6.9) 


Substituting equation (6.8) into (6.9) and using q as generalized coordinates the actuator 
equation becomes 


where B comes from equation (5.29) and Ca is the applied voltage to the piezoceramic 
element. This is the modified version of equation (6.3) that must be used for the 
dynamical system model. 

If the voltage e is taken in terms of generalized coordinates, the equation in terms 
of voltage is 


( 6 . 11 ) 


Equation (6.11) represents the voltage output from a piezoceramic sensor. 
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Using the material properties from Table 6.1 and the preceding equations of 
motion, the finite element model produced twelve natural modes of vibration. Appendix 
A contains the MATLAB script file, fem.m, that was used to solve this model. Table 6.2 
gives the values for the natural frequencies and Figures 6.5 and 6.6 give the modeshapes 
for the first two modes. Appendix B contains the graphs for all twelve modeshapes. 
These two modes are the primary carriers of energy for the structure and it is they that the 
compensator will be designed to actively control. 


Mode 


1 

0. 2958 

2 

0.8707 

i 3 

11.108 

4 

28.496 

5 

45.144 

i_ 6 _ 

102.78 

! 7 

278.83 

8 

341.04 

9 

1347.7 

10 

1392.7 

11 

5284.5 

12 

5343.0 


Table 6.2. Finite Element Model 
Natural Frequencies 


Convergence of a finite element model to its analytical counterpart is dependent 
on the number of elements used in the model, the higher the number of elements the 
quicker the model is to converge to the actual values. A finite element model’s accuracy 
can be roughly estimated as half the number of elements used. In this case, six elements 
are used and thus the first three modes should be fairly accurate. Using elementary modal 
analysis, the first three modes of vibration are known to be 0.288 Hz , 0.877 Hz, and 9.82 
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Hz, respectively. This indicates percentage errors for the first three modes as 2.7%, 
0.72%, and 13.1 % each. The diminishing accuracy with the increase in mode number is 
readily apparent due to a small number of elements in the model. 



Figure 6.5. Modeshape for First Mode 



At this time it should be noted that the first finite element model consisted of 
twelve elements giving twenty-four natural modes. This first model’s accuracy extended 
to the sixth natural mode before significant divergence was noticed. The problem 
encountered with the higher fidelity model was the introduction of uncontrollable and 
unobservable states. As discussed in Chapter m, the solution was to simplify the model 
to the current twelve mode model for the remainder of the thesis work. 

Since the piezoceramic actuators use electro-mechanical coupling from an applied 
voltage to a strain on the beam, only a moment, or torque, may be applied as the control 
input. Based on the given modeshapes (eigenvectors) and the knowledge that the region 
of greatest torque is at the base for a fixed-free (cantilevered) configuration, it can be 
inferred that the optimal actuator and sensor placement be at the base of the structure and 
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at the base of the second arm. This configuration will give the greatest control input 
possible and also provide a nearly collocated sensor/actuator pair. 

C. TRANSFORMATION TO MODAL FORM 

Equation (6.10) is a damped free equation of motion derived from finite element 
analysis incorporating piezoceramic electro-mechanical coupling. Actually the FSS 
flexible arm has natural structural damping and in order to increase the fidelity of the 
model of our structure, the natural damping must be included. In order to add in the 
structural damping a linear transformation of the state must be used to get the model into 
modal form. Modal form is one that de-couples the equations into linearly independent 
first-order differential equations. In Chapter 111, Section C., a method for a linear 
transformation of state is given. Using the generalized eigenvalue approach, there is a 
transformation matrix, O, that transforms equation (6.10) into diagonal form. The 
response of the sy.stem, in the absence of an external input is given by 

Mq + Kq = 0 (6.12) 

The transformation from physical or generalized coordinates, q, to modal 
coordinates, is detennined by the linear similarity transformation 

q = cl)^- T = <£-' ^ = Tq (6.13) 

The desired equations of motion are of the form 
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Mq + Cq + Kq = 0 


(6.14) 


where C is the damping matrix for the system in physical coordinates. Combining 
equation (6.13) with equation (6.14) we obtain 


M®i + Cd)^ + K<I><^ = 0 

+ cD’'K<D<J = 0 


(6.15) 


where F is determined so that 


(fMO = I 


(D’^CO = 


2 ^ 6 ), 






In this way, the equations of motion in modal coordinates become a set of de - coupled 


second - order differential equations 


^ + [fi]^ + [A]| = 0 (6.16) 

which can be written in state - space form as 


^ = -[Q]^-[A]^ 



Using data obtained from previous experimentation on the FSS flexible arm, it 
was determined that the first two natural modes need 0.4% as the damping ratio and 10% 


46 



for the higher modes. These values are inserted into equation (6.17) and the equation is 
then transformed back to physical coordinates. Recalling equation (6.13c) a state-space 
transformation matrix can be obtained 


= [T]q 


(6.18) 


thus the system matrix transformed back to physical coordinates is 

[T]ij = A„Tq 

q = [T]'A„[T]q ^ q = Aq (6.19) 

where A=[T]''a„[T] 


D. STATE-SPACE REPRESENTATION 

The analytical model of our system can be represented by the following state- 
space first-order coupled equations 


it = Ax + Bu + Fw 
y = Cx + V 


(6.20) 


The system matrix. A, from equation (6.19) is now an accurate twenty-four state 
model approximating the actual FSS flexible arm assembly. The control matrix, B, based 
on equations (5.29) and (6.10), has two independent actuators and is a 24x2 matrix. The 


47 




sensor matrix, C, is a 6x24 matrix that contains two piezoceramic sensor outputs and four 
outputs from an optical infrared VisionServer camera system developed by Real time 
Innovations, Inc. [Ref. 6] The piezoceramic outputs are governed by equations (5.29) and 
(6.11). The VisionServer measures the displacements and rotations of the elbow and tip 
assemblies individually and outputs those four states directly. The plant uncertainty 
matrix, F, is initially the identity matrix and the model uncertainty is approximated at 5%. 
The sensor noise is taken as the squares of component ims noise values from 
documentation and previous research 




Vn. ANALYTICAL MODEL SIMULATION 


In this chapter a LQG controller is designed using the methods discussed in 
Chapter El. The LQG controller will be designed to meet given performance criteria 
outlined in Section B of this chapter. The analytical model will be simulated using 
MATLAB and SIMULUMK after which it will be ported to a VAXstation running 
MATRIXx to be tested on the real-time controller. 

A. MODEL REPRESENTATION 

From equation (6.20) the model is represented in state-space form as 

X = Ax + Bu + Fw 
y = Cx + v 


The model is formed by the following components: the flexible body dynamics, 
the piezoceramic actuators, the piezoceramic sensors, and the VisionServer camera 
system[Ref. 6]. The type of controller to be designed is a regulator, i.e., there is no 
reference or commanded input to track, since the objective is to dampen all the 
oscillations in the system. The major disturbances are in the form of continuous and 
impact vibrations to the flexible arm, while secondary disturbances are the noise inputs to 
the sensors and actuators from the control system hardware. The state variables 
considered in the nominal model are the nodal displacements and nodal rotations, w and 
0. The model upon which the LQG is designed combines the nominal plant, the 
estimated plant uncertainties, and the estimated process-noise model. 
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B. LQG DESIGN REQUIREMENTS 


The primary objective is to damp out the effect of the disturbances as quickly as 
possible and maintain stability robustness. Specifically, the most important variable to 
minimize is the rotational displacement of the tip of the arm since this directly affects the 
performance of the antenna structure in space. 

The FSS is setup for an initial condition response to simulate an impact 
disturbance input. The piezoceramic actuators are low authority controllers, and the time 
constant for settling out the disturbance will be on the order of five seconds. The steady 
state error to a continuous sinusoidal input is to be less than 0.1 ° for the tip assembly. In 
summary the performance objectives for the LQG controller are 

Settling Time Tj < 20 seconds (Time constant = 5 seconds) 

Steady State Error e* < 0.1° (Due to a sinusoidal input) 


C. PROPOSED DESIGN 

The design of the LQG controller for this model is given by a combination of the 
solutions of the Linear Quadratic Regulator (LQR) and the Kalman Estimator (LQE) 
problem. The optimal control solution of the LQR problem which minimizes the 
performance index given by equation (3.6) is given by the state feedback equation 


u = -K 




(7.1) 
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where Klqr is given by equation (3.7) provided that the state weighting matrix, Q, is 
symmetric and positive semi-definite, the control weighting matrix, R, is symmetric and 
positive definite, and the pair (A,B) is controllable, with the controllability matrix given 
by equation (3.14a). Since the sensors do not measure all the states, it is desirable to 
transform equation (3.6) into an output dependent equation so as to have direct influence 
over the weighting parameters. The cost function, in general given as 

J= j(x’^Qx+u’^Ru)rfr 


can be adapted to the case when we penalize the output, y. Using the fact that y=Cx and 
y'^=x'^C^ , another matrix, Q^can be defined where C’^Q^C, and substitution into the 
above equation yields 


•/= J(y"Q>,y+u^R‘*)‘^/ (7-2) 

which has an explicit dependence on the sensor outputs. 

The solution to the LQR seeks a compromise between minimum energy (control 
input) and best performance. For an antenna structure in space it is desirable to minimize 
the beam losses incurred due to the disturbance inputs. Similarly, in this model the 
priority is to minimize the tip deflection and rotation so we will penalize those outputs 
more heavily than the others. In order to keep the output of the piezoceramic sensors 
within limitations a penalty must be imposed on those outputs as well. 
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Using SIMULINK, a dynamical model was created and simulated. Figure (7.1) is 
the block diagram of the simulation model. 



Figure 7.1. SIMULINK Block Diagram 


D. SIMULATION AND RESULTS 

From simulation, the values of Qy, the output weighting matrix, and R, the control 
weighting matrix, were iteratively determined to be 

100 0 0 0 0 O' 

0 100 0 0 0 0 

0 0 1 0 0 0 ro.1 Ol (7.3) 

^'" 0 0 0 1 0 0 **"[o 0.lj 

0 0 0 0 9000 0 

0 0 0 0 0 900000 

where the sensor output vector, y, contains the piezoceramic sensors as the first two 
elements and four VisionServer outputs (elbow displacement, elbow rotation, tip 
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displacement, and tip rotation) as the last four elements of the output vector. The two 
control inputs are the base actuator and the elbow actuator, respectively. 

These values kept the control inputs within their limitations of +150 volts, the 
sensors within their limitations of ±10 volts, minimized the steady state error, and met the 
settling time constraints. Figure 7.2 illustrates a typical tip rotational displacement in the 
time domain, it is representative of all sensor data observed. Appendix C contains the 
simulation results for all the sensors. 

Another method for determining the performance of the LQG controller is to 
measure the damping effectiveness for the system. A good method for determining the 
damping from simulations is the Log Decrement method given by 


? = :rrin(^) 


(7.4) 


with ^ being the damping ratio, n the number of cycles between measurements, and At 
the initial and final amplitudes, respectively. In order to determine damping, we need to 
excite one mode at a time. In general, we were able to excite independently only the first 
two modes. 
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Open Loop Response to Initial Condition 



a.-1,5‘-—'-'-'-'-'-' 

F 0 5 10 15 20 25 30 


LQR Response to Initial Condition 



time (seconds) 


Figure 7.2. Tip Rotational Displacement Initial Condition Response 

The damping coefficient for the first two modes determined from the simulation 
results by the log decrement method is 12%. Recalling the original damping ratio of 
0.4% an increase in damping on the order of 29 times or 2,900% is obtained. 

The next step is to design the Kalman filter. The solution to the optimal state 
estimator is given by equation (3.9) where the observer gain is taken from equation 
(3.11). Equation (3.11) has a solution provided that W, the plant uncertainty variance 
matrix, is symmetric and positive semidefinite; the sensor noise variance matrix, V, is 
symmetric and positive definite, and the pair (A,C) is observable where the observability 
matrix is given by equation (3.14b). The values determined from hardware 
documentation for the sensors, and previous research for the plant uncertainty model and 
process-noise model are given below. Of note is the fact that the piezoceramic sensors 
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generate a clean signal that has very little noise and that is reflected in the corresponding 
variance matrix. 


w = o.oi[:i:] 


"O.OOOl 0 0 0 0 0 " 

0 0.0001 0 0 0 0 

0 0 0.00025 0 0 0 

0 0 0 0-002 0 0 

0 0 0 0 0.00025 0 

0 0 0 0 0 0.001_ 


(7.4) 


The designs for the LQR controller and the Kalman filter were been performed 
individually. To create the compensator for the system, we must combine the two 
together into one system. Using equation (3.9), the estimator equation, we will introduce 
the controller equation (7.1) for the control input modified for the estimated states rather 
than the actual states 


X = Ax + Bu + L(y - Cx) 


Substituting the control equation into the estimator equation we arrive at the 
overall compensator equation 


i = (A-BK,Q„-LC)x + Ly 


(7.5) 
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with the inputs to the compensator being the sensor outputs, y, and the outputs of the 
compensator as the control input to the structural system ,u. 

The response using this LQG controller can be expected to converge to the LQR 
response over time. The speed of convergence is a reflection of the estimator’s time 
constant. Figure 7.3 illustrates the error between the LQR response and the LQG 
response due to the Kalman filter. From Figure 7.3 it can be seen that the LQG response 
is less effective than the LQR due to the loss of robustness with the introduction of the 
uncertainty and noise into the system. 



Figure 7.3. Convergence of LQG to the LQR Initial Condition Response. 

The value for the initial condition, the tip displaced by 3 centimeters, was 
determined experimentally by defining the maximum deflection that would keep the 
voltage input to the controller by the piezoceramic sensors within limits (±10 volts). 

The final step in the design process for the LQG controller is to discretize the 
system for the experimental verification effort. The method of discretization to be used is 
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the Tustin transformation. The relationship between the continuous domain (s-domain) 
and the discrete domain (z-domain) is given by 

z = (7.6) 

where T is the sampling period. A Tustin transformation maintains the frequency 
response of the continuous system while preserving the mapping of the s-plane into the z- 
plane The Tustin transformation is a Fade approximation to the exponential, equation 
(7.6). [Ref 7, pp. 253-282] 


sT 

1 + — 



The continuous system will be transformed to the discrete domain using a sampling 
period of 0.01 seconds (100 Hz). The tip response for the LQG controller in the discrete 
domain is shown in Figure 7.4, note the similarity between the continuous response and 
the discrete response. 

The implementation on the real hardware is shown in [Ref 8]. The controller is 
implemented in software on the VAXstation 3100 using MATRIXx and SystemBuild. 
The real-time code is then downloaded to the AC-100 real-time controller for 
experimental verification of this design. 
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LQG 



Figure 7.4. Comparison Between the LQG Continuous and Discrete Responses 
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VIII. CONCLUSIONS AND RECOMMENDATIONS 


A. SUMMARY 

This thesis presents the design of an optimal controller for a flexible spacecraft 
structure model using Linear Quadratic Gaussian techniques with piezoceramic actuators 
and sensors. An optical sensor package is also introduced into the system to provide 
control inputs. 

A finite element model of the structure was provided that accounted for the 
electro-mechanical coupling effects of the piezoceramic devices. The finite element 
model was enhanced by the addition of structural damping. The design of the Linear 
Quadratic Gaussian controller is based on the combination of the Linear Quadratic 
Regulator (LQR) solution and the optimal estimator problem (Kalman filter). The 
solution to the LQR problem was transformed from a state dependent problem to a output 
dependent problem. This technique allowed for the direct control over the control inputs 
in relation to the sensor outputs. The plant uncertainty and process-noise matrices were 
determined by iteration using the sensor specifications and previous research as the 
starting points and then fine-tuning each of the matrices’ elements to achieve the desired 
performance. The estimator and LQR controller were combined to form the LQG 
controller and subjected to rigorous continuous-time simulation in Matlab and 
SIMULINK. 

The model was then transformed to discrete space via a Tustin transformation 
following the continuous-time simulations. The discrete controller was also simulated in 
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Matlab to ensure a stable introduction into the experimental setup with the real-time 
controller, AC-100. 


B. CONCLUSIONS 

In this study a Linear Quadratic Gaussian controller was designed and developed 
for the Flexible Spacecraft Simulator (FSS) at the Naval Postgraduate School. All 
specifications for the FSS were met to include the sensor and actuator limitations for the 
real-time controller. The controller responded to a initial condition disturbance and met 
the time requirements and steady state eiror requirements specified in the design criteria. 
A significant increase in damping on the order of 2900% was achieved with the LQG 
controller. The optimal estimator’s performance met specifications, but only marginally 
due to the fact that some of the states had small observability values. This condition may 
have severe control implications when implemented experimentally. 

C. RECOMMENDATIONS 

The LQG controller design should be verified experimentally by testing with the 
FSS flexible structure [Ref. 8]. The finite element method for characterizing the structure 
is good for approximation only, and the use of modal analysis would increase the fidelity 
of the model. Ftuther study in the use of Loop Transfer Recovery techniques for this 
controller can recover some of the loss in robustness, but this is at the price of losing 
nominal perfonnance. 
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APPENDIX A 


FINITE ELEMENT MODEL MATLAB CODE 

%%Finite EleMent Model PrograM for the FSS 
%% 

%%This prograM characterizes the FSS flexible appendage. 

%%LCDR Bill Harrington, USN October 1995 

%%The flexible appendage is to be Modeled as a 6 eleMent FEM with 
%%piezoceraMic actuators located on eleMents 1 and 4, and piezoceraMic 
%%sensors located on eleMents 2 and 5. The elbow and tip will be treated ; 
%%rigid body Masses. 


%%This Model assuMes 0.5% structural daMping, and will solve for the first 
%%] 2 natural Modes. The degree of accuracy is suitable only for the 
%%first 3 Modes (half the nuMber of eleMents, normally). 


clc 

format short e 
global A B C D 

%%stnictural properties 

nuMeleMents=6; 

h=[0.07 0.03 0.5688 0.07 0.03 0.549]; 
tlncKness=0.0015875; 
heighl=0.0254; 
rho=2800; 

elMass=rho*height*thicK,ness: 
el inertia=l /12 * thicKness^3 ^height; 
E=72e9; 


%%Meters 


%%Kg/M''3 

%%Mass/lengtlr 

%%Meter''4 

%%Modulus of elasticity 


%%CoMpute eleMental stif&iess and Mass Matrices 


KeleMent=0*ones(4,4*nuMeleMents); 

MeleMent=0*ones(4,4*nuMeleMents); 


for i=l:nuMeleMents,... 

KeleMent(:,4*i-3:4*i)=E*elinertia/h(i)''3*[12 6*h(i) -12 6*h(i);... 
6*h(i) 4*h{i)''2 -6*h(i) 2*h(i)^2;-12 -6*h(i) 12 -6*h(i);... 

6*h(i) 2*h(i)^2 -6*h(i) 4*h(i)''2]; 


for i-^1 muMeleMents,... 

MeleMent(:,4*i-3:4*i)=elMass*h(i)/420*[156 22*h(i) 54 -]3»h(i);... 
22*h(i) 4*h(i)^2 13*h(i)-3*h(i)''2;54 13*h(i) 156 -22*h(i);... 
-13*h(i) -3*h(i)''2 -22*h(i) 4*h(i)''2]; 
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%%Construct global Mass and stif&ess Matrices 

M=0*ones(2*(nuMeleMents+l)); 

K=0*ones(2*(nuMeleMents+l)); 

for i=l inuMeleMents,... 

M(2*i-l:2*i+2,2*i-l:2*i+2)=M(2*i-l;2*i+2,2*i-l:2*i+2)+MeleMent(:,4*i-3;4»i); 

end 

for i=l:nuMeleMents,... 

K(2*i-1 ;2*i+2,2*i-1 ;2*i+2)=K(2*i-l :2*i+2,2*i-1:2*i+2) +KeleMent(:,4*i-3:4»i); 
end 

%%Clean up 
clear elinertia; 

%%Fix for fixed boundry condition (row,coluMns 1&2 =0) 

M=M(3:2*(nuMeleMents+1 ),3;2*(nuMeleMents+1)); 

K=K(3 ;2*(nuMeleMents+l ),3:2*(nuMeleMents+1)); 

M(5,5:8)=M(5,5:8)-MeleMent( 1,13:16); 

M(6:8,5)=M(6;8,5)-MeleMent(2:4,13); 

K(5,5:8)=K(5,5:8)-KeleMent( 1,13:16); 

K(6:8,5)=K(6:8,5)-KeleMent(2:4,13); 


clear KeleMent MeleMent; 

%%Now add elbow and tip Masses, 

Melbow=0.40823; 

Mtip=0,37648; 

MbeaM2=(h( 1 )+h(2)+h(3))*elMass; 

M(5,5)=M(5,5)+Melbow+MbeaM2+Mtip; 

M(11,11)=M(11,11)+Mtip; 

clear Mpoint Melbow Mtip ipoint MbeaM2; 

%%Define the piezo actuator and sensor eleMents 

delta=thicKness/2; 

tp=2*1.905e-4; 

wp=0.02; 

Ep=6.3el0; 

d31=-1.8e-10; 

eT3=1.5e-8; 

rhop=7700; 

Kpiezo=wp*tp*Ep*(delta^2 + delta*tp + (tp"2)/3); 

Ka=Kpiezo/h(l)»[12/h(l)''2 6/h(l) -12/h(l)''2 6/h(l);... 
6/h(l) 4 -6/h(l) 2;-12/h(l)''2 -6/h(l) 12/h(l)^2 -6/h(l);... 
6/h(l)2-6/h(l)4]; 
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Ks=Kpiezo/h(2)*[12/h(2)''2 6/h(2) -12/h(2)''2 6/h(2);... 

6/h(2) 4 -6/h(2) 2;-12/h(2)''2 -6/h(2) 12/h(2)''2 -6/h(2);... 

6/h(2) 2 -6/h(2) 4]; 

Mp=rhop*wp*tp; %%Mass/length piezo 

Ma=Mp»h(l)/420*[156 22*h(l) 54 -13*h(l);22*h(l) 4*h(ir2 13*h(l) -3*h(ir2;... 
54 13*h(l) 156 -22*h(l);-13*h(l) -3*h(ir2 -22»h(l) 4*h(l)''2]; 

Ms=Mp*h(2)/420*[156 22»h(2) 54 -13*h(2);22*h(2) 4*h(2r2 I3*h(2) -3»h(2)''2;... 
54 13»h(2) 156 -22*h(2):-13*h(2) -3*h(2r2 -22*h(2) 4*h(2r2]; 

%%Add the piezo eleMent Mass and stiffiiess Matrices to the structure 

M( 1:2,1:2)=M( 1:2,1:2)+2*Ma(3:4,3:4); 

M( 1 ;4,1:4)=M( 1:4,1:4)+2*Ms; 

M(6:8,6:8)=M(6:8,6:8)+2*Ma(2:4,2:4); 

M(7:10,7:10)=M(7:10,7:10)+2*Ms; 

K(1:2,1:2)=K(1:2,1:2)+2*Ka(3:4,3:4); 

K(l:4,l:4)=K(l:4,l:4)+2*Ks; 

K(6:8,6:8)=K(6:8.6:8)+2*Ka{2:4,2:4); 

K(7:10,7:10)=K(7:10,7:10)+2*Ks; 

clear Ks Ka Ms Ma Kpiezo; 

%%Solve for natural frequencies and Mode shapes. 

[oMega2,Phi,Psi]=eign(K,M); 

oMega=sqrt(oMega2); 

Hertz-oMega/2/pi; 

ttl=str2mat(' Oniega','Hertz',' '); 

[ttl(l,:)ttl(3,:)ttl(2,:)] 

[oMega Hertz] 

%%Construct the a,b,c, and d Matrices for state-space forM. 

%%Need More info on the piezos. 


bbl=0; bb2=-d31»Ep»wp*(delta +tp/2); bb3=0; bb4=-bb2; 
gaMMa=wp»h(2)/tp*(eT3-d31''2*Ep); 
bl=[bb3bb4 0 000000000]; 
b2=[0 0 0 0 bbl bb2 bb3 bb4 0 0 0 0]; 

btMp=[bl' b2']; btMp2=-2»inv(M)*btMp; btMp3=-2*Phi'*btMp; 
ctMp=[bbl bb2 bb3 bb4 zeros( 1,20)]/gaMMa; %sensor 1 

ctMp2=[0 0 0 0 0 0 bb 1 bb2 bb3 bb4 0 0 zeros] 1,12)]/gaMMa; “/osensor 2 

ctMp3=[0 0 0 0 1 0 zeros] 1,18)]; %e. disp. 

ctMp4=[0 0 0 0 0-1 zeros] 1,18)]; %e. rot. 

ctMp5=[zeros] 1,10) 1 0 zeros] 1,12)]; %t. disp. 

ctMp6=[zeros] 1,10) 0 -1 zeros] 1,12)]; %t. rot. 


a=[0*ones]12) eye]12);-inv]M)*K 0*ones]12)]; 
b=t0*ones] 12,2);btMp2]; 

c=[ctMp;ctMp2;ctMp3;ctMp4;ctMp5;ctMp6]; %observation matrix 

d=0*ones]6,2); 
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clear bbl bb2 bb3 bb4 gaMMa bl b2 btMp ctMp; 
clear ctMp2 ctMp3 ctMp4 ctMpS ctMp6 ctMp7 ctMp8; 
clear E Ep eT3 i; 

clear btMp2 cx d31 delta e elMass ep et3 height Mp rho rhop; 
clear thicKness tp wp 

%%Construct transformation matrix from modal to physical coordinates. 

Phi_inv=inv(Phi); 

PPhi=zeros(24); 

PPhi(l:12,l:12)=PhiJnv; 

PPhi( 13:24,13:24)=Phi_inv; 

PPhiJnv=inv(PPhi); 

%%Insert natural structural damping, 0.5% for first two modes then 20% for rest. 

damp=-2*0.005*oMega( 1:2);damp2=-2*0.1 *oMega(3:12); damp3=diag([damp' damp2']); 

a2=[zeros(12) eye( 12);-diag(oMega2) damp3]; 

b2=[zeros(12,2);btMp3]; 

c2=[c(:,l:12) zeros(6,12)]*[Phi zeros(size(Phi));zeros(size(Phi)) Phi]; 
d2=d; 


a2m=PPhi_inv*a2*PPhi; 

A=clean(a2m,le-6); 

b2m=PPhi_inv*b2; 

B=clean(b2m, 1 e-6); 
c2m=c2*PPhi; 

C=clean(c2m, 1 e-7); 

D=d; 

%%Set an initial condition (tip displaced by 3 cm) 
x0=[0 0 0 0 0.02 0 0 0 0 0 -0.03 0 zeros(l,12)]; 

%%Time vector, t 
t=0;0.01:30; 

dispCInitialization of variables complete for two beam analysis. (FEM3)') 
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APPENDIX B 


ANALYTICAL MODEL MODESHAPES 
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APPENDIX C 


SENSOR OUTPUTS FOR LQG CONTROLLER 
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LQG Response to Initial Condition 
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